Source-transformation Differentiation of a C++-like Ice Sheet model

International audienceAlgorithmic Differentiation (AD) has become one of the most powerful tools to improve our understanding of theEarth System. If AD has been used by the ocean and atmospheric circulation modelingcommunity for almost 20 years, it is relatively new in the ice sheet modeling community. The Ice SheetSystem Model (ISSM) is a C++, object-oriented, massively parallelized, new generation ice sheet model that recentlyimplemented AD to improve its data assimilation capabilities. ISSM currently relies on Object Overloading throughADOL-C and AMPI. However, experience shows that Object Overloading AD on ISSM is significantly more memoryintensive compared to the primal code. We want to investigate other AD approaches to improve the performance ofthe AD adjoint. Yet, to our knowledge, there is no source-to-source AD tool that supports C++.To overcome this problem, we have developed a prototype of ISSM entirely in C, called Boreas, in order to testsource-to-source transformation and compare the performance of these two approaches to AD. Boreas is a clone ofISSM, the main difference with ISSM is that all the objects are converted to C-structures and some function nameshave been adapted in order to be unique, but the code architectures are identical. The programming style of Boreas isa first attempt at defining a programming style of (or a sub-language of) C++ that source-transformation AD couldhandle. We present here the first results of Source-Transformation AD of Boread with the AD tool Tapenade

Extended enthalpy formulations in the Ice-sheet and Sea-level System Model (ISSM) version 4.17: discontinuous conductivity and anisotropic streamline upwind Petrov--Galerkin (SUPG) method

The thermal state of an ice sheet is an important control on its past and future evolution. Some parts of the ice sheet may be polythermal, leading to discontinuous properties at the cold–temperate transition surface (CTS). These discontinuities require a careful treatment in ice sheet models (ISMs). Additionally, the highly anisotropic geometry of the 3D elements in ice sheet modelling poses a problem for stabilization approaches in advection-dominated problems. Here, we present extended enthalpy formulations within the finite-element Ice-Sheet and Sea-Level System model (ISSM) that show a better performance than earlier implementations. In a first polythermal-slab experiment, we found that the treatment of the discontinuous conductivities at the CTS with a geometric mean produces more accurate results compared to the arithmetic or harmonic mean. This improvement is particularly efficient when applied to coarse vertical resolutions. In a second ice dome experiment, we find that the numerical solution is sensitive to the choice of stabilization parameters in the well-established streamline upwind Petrov–Galerkin (SUPG) method. As standard literature values for the SUPG stabilization parameter do not account for the highly anisotropic geometry of the 3D elements in ice sheet modelling, we propose a novel anisotropic SUPG (ASUPG) formulation. This formulation circumvents the problem of high aspect ratio by treating the horizontal and vertical directions separately in the stabilization coefficients. The ASUPG method provides accurate results for the thermodynamic equation on geometries with very small aspect ratios like ice sheets

Greenland Subglacial Discharge as a Driver of Hotspots of Increasing Coastal Chlorophyll Since the Early 2000s

Subglacial discharge emerging from the base of Greenland\u27s marine-terminating glaciers drives upwelling of nutrient-rich bottom waters to the euphotic zone, which can fuel nitrate-limited phytoplankton growth. Here, we use buoyant plume theory to quantify this subglacial discharge-driven nutrient supply on a pan-Greenland scale. The modeled nitrate fluxes were concentrated in a few critical systems, with half of the total modeled nitrate flux anomaly occurring at just 14% of marine-terminating glaciers. Increasing subglacial discharge fluxes results in elevated nitrate fluxes, with the largest flux occurring at Jakobshavn Isbræ in Disko Bay, where subglacial discharge is largest. Subglacial discharge and nitrate flux anomaly also account for significant temporal variability in summer satellite chlorophyll a (Chl) within 50 km of Greenland\u27s coast, particularly in some regions in central west and northwest Greenland

A JavaScript API for the Ice Sheet System Model: towards on online interactive model for the Cryosphere Community

Abstract. Earth System Models (ESMs) are becoming increasingly complex, requiring extensive knowledge and experience to deploy and use in an efficient manner. They run on high-performance architectures that are significantly different from the everyday environments that scientists use to pre and post-process results (i.e. MATLAB, Python). This results in models that are hard to use for non specialists, and that are increasingly specific in their application. It also makes them relatively inaccessible to the wider science community, not to mention to the general public. Here, we present a new software/model paradigm that attempts to bridge the gap between the science community and the complexity of ESMs, by developing a new JavaScript Application Program Interface (API) for the Ice Sheet System Model (ISSM). The aforementioned API allows Cryosphere Scientists to run ISSM on the client-side of a webpage, within the JavaScript environment. When combined with a Web server running ISSM (using a Python API), it enables the serving of ISSM computations in an easy and straightforward way. The deep integration and similarities between all the APIs in ISSM (MATLAB, Python, and now JavaScript) significantly shortens and simplifies the turnaround of state-of-the-art science runs and their use by the larger community. We demonstrate our approach via a new Virtual Earth System Laboratory (VESL) Web site

Drivers of Change of Thwaites Glacier, West Antarctica, Between 1995 and 2015

We run several transient numerical simulations applying these three perturbations individually. Our results show that ocean-induced ice-shelf thinning generates most of the observed grounding line retreat, inland speed-up, and mass loss, in agreement with previous work. We improve the agreement with observed inland speed-up and thinning by prescribing changes in ice-shelf geometry and a reduction in basal traction over areas that became ungrounded since 1995, suggesting that shelf breakups and thinning-induced reduction in basal traction play a critical role on Thwaites's dynamics, as pointed out by previous studies. These findings suggest that modeling Thwaites's future requires reliable ocean-induced melt estimates in models that respond accurately to downstream perturbations

Implementation and performance of adaptive mesh refinement in the Ice Sheet System Model (ISSM v4.14)

Accurate projections of the evolution of ice sheets in a changing climate require a fine mesh/grid resolution in ice sheet models to correctly capture fundamental physical processes, such as the evolution of the grounding line, the region where grounded ice starts to float. The evolution of the grounding line indeed plays a major role in ice sheet dynamics, as it is a fundamental control on marine ice sheet stability. Numerical modeling of a grounding line requires significant computational resources since the accuracy of its position depends on grid or mesh resolution. A technique that improves accuracy with reduced computational cost is the adaptive mesh refinement (AMR) approach. We present here the implementation of the AMR technique in the finite element Ice Sheet System Model (ISSM) to simulate grounding line dynamics under two different benchmarks: MISMIP3d and MISMIP+. We test different refinement criteria: (a) distance around the grounding line, (b) a posteriori error estimator, the Zienkiewicz-Zhu (ZZ) error estimator, and (c) different combinations of (a) and (b). In both benchmarks, the ZZ error estimator presents high values around the grounding line. In the MISMIP+ setup, this estimator also presents high values in the grounded part of the ice sheet, following the complex shape of the bedrock geometry. The ZZ estimator helps guide the refinement procedure such that AMR performance is improved. Our results show that computational time with AMR depends on the required accuracy, but in all cases, it is significantly shorter than for uniformly refined meshes. We conclude that AMR without an associated error estimator should be avoided, especially for real glaciers that have a complex bed geometry.121215232CONSELHO NACIONAL DE DESENVOLVIMENTO CIENTÍFICO E TECNOLÓGICO - CNPQ140186/2015-

The transferability of adjoint inversion products between different ice flow models

Among the most important challenges faced by ice flow models is how to represent basal and rheological conditions, which are challenging to obtain from direct observations. A common practice is to use numerical inversions to calculate estimates for the unknown properties, but there are many possible methods and not one standardised approach. As such, every ice flow model has a unique initialisation procedure. Here we compare the outputs of inversions from three different ice flow models, each employing a variant of adjoint-based optimisation to calculate basal sliding coefficients and flow rate factors using the same observed surface velocities and ice thickness distribution. The region we focus on is the Amundsen Sea Embayment in West Antarctica, the subject of much investigation due to rapid changes in the area over recent decades. We find that our inversions produce similar distributions of basal sliding across all models, despite using different techniques, implying that the methods used are highly robust and represent the physical equations without much influence by individual model behaviours. Transferring the products of inversions between models results in time-dependent simulations displaying variability on the order of or lower than existing model intercomparisons. Focusing on contributions to sea level, the highest variability we find in simulations run in the same model with different inversion products is 32 , over a 40-year period, a difference of 3.67 mm. There is potential for this to be improved with further standardisation of modelling processes, and the lowest variability within a single model is 13 , or 1.82 mm over 40 years. While the successful transfer of inversion outputs from one model to another requires some extra effort and technical knowledge of the particular models involved, it is certainly possible and could indeed be useful for future intercomparison projects